Implementation of the Jacobian-free Newton-Krylov

14 We have implemented the Jacobian-free Newton-Krylov (JFNK) method 15 for solving the first-order ice sheet momentum equation in order to improve 16 the numerical performance of the Community Ice Sheet Model (CISM), the 17 land ice component of the Community Earth System Model (CESM). Our 18 JFNK implementation is based on significant re-use of existing code. For ex19 ample, our physics-based preconditioner uses the original Picard linear solver 20 in CISM. For several test cases spanning a range of geometries and bound21 ary conditions, our JFNK implementation is 1.84-3.62 times more efficient 22 than the standard Picard solver in CISM. Importantly, this computational 23 gain of JFNK over the Picard solver increases when refining the grid. Global 24 convergence of the JFNK solver has been significantly improved by rescaling 25 the equation for the basal boundary condition and through the use of an in26 exact Newton method. While a diverse set of test cases show that our JFNK 27 implementation is usually robust, for some problems it may fail to converge 28 with increasing resolution (as does the Picard solver). Globalization through 29 parameter continuation did not remedy this problem and future work to im30 prove robustness will explore a combination of Picard and JFNK and the use 31 of homotopy methods. 32